Burgers' equation¶
%matplotlib inline
%config InlineBackend.figure_format = 'svg'
from utils.snapshot_widgets import interact
from ipywidgets import widgets
from ipywidgets import FloatSlider
from exact_solvers import burgers
from exact_solvers import burgers_demos
In this chapter, we study another scalar nonlinear conservation law: Burgers' equation. Burgers' equation models momentum transport in a fluid of uniform density and pressure, and it is the simplest equation that captures some key features of gas dynamics.
To examine the Python code for this chapter, see:
Burgers' equation has been used extensively for developing both theory and numerical methods, and it will allow us to further explore the Riemann problem for nonlinear conservation laws. The equation is given by
\begin{align}
q_t + \left(\frac{1}{2}q^2\right)_x = 0.
\label{burgers0}
\end{align}
It is a scalar conservation law with the nonlinear flux term $f(q)=q^2/2$. The quasilinear form is obtained by applying the chain rule to differentiate the flux term:
\begin{align*}
q_t + qq_x = 0.
\end{align*}
This equation looks very similar to the advection equation, with the difference that the advection speed at each point is given by the solution $q$. The nonlinearity in the Burgers flux term is essentially the same as in the momentum conservation equation in fluid dynamics, so studying its dynamics provides a guideline to understand more complex fluid dynamics, as we will do in Euler.html.
Shock formation¶
As the speed of propagation in Burgers' equation is given by $q$ itself, the peak of a wave travels faster than the rest. This behavior is illustrated in the next figure. The dashed line shows the initial condition while the solid lines show the solution at later times.
interact(burgers_demos.multivalued_solution,
t=FloatSlider(min=0.,max=9.,value=7.75));
Notice that at first $q$ remains single-valued for every $x$. However, after some time the crest of the wave overtakes the leading edge. After this time, we obtain a triple-valued solution for certain values of $x$, which is unphysical since momentum is single-valued. The first time this overtaking happens is referred to as the breaking time -- a reference to waves breaking on the beach. It is also the point where the conservation law, in its differential form, breaks down and where the characteristics cross each other for the first time. As we learned in Traffic_flow.html, when characteristics cross, a shock wave forms. Mathematically, imposing a shock will avoid the problem of the solution being multivalued. Where should the shock be located?
If we replace part of the multivalued solution interval with a shock, some mass will be removed (area $A_1$ in the figure below) and some mass will be added (area $A_2$ in the figure below). In the live notebook, the figure below shows possible locations for the shock at a given time; in order to maintain conservation of the integral of $q$, the shock must be placed so that these two areas are equal.
interact(burgers_demos.shock_location,
xshock = FloatSlider(min=6,max=10,step=0.25,value=7.75));
This geometric reasoning provides a nice intuition for the shock location, but is cumbersome in practice. To determine the path of the shock we can again use the Rankine-Hugoniot jump condition, derived in Traffic_flow.html. Since the flux for Burgers' equation is $f(q) = q^2/2$, we have
\begin{align*}
s(q_r-q_l) & = f(q_r) - f(q_l)\\
& = \frac{1}{2} (q_r^2 - q_l^2) \\
\Rightarrow \ \ \ \ s &= \frac{1}{2}(q_l + q_r).
\end{align*}
Here $q_l, q_r$ denote the solution values just to the left and right of the shock. In general (as in the image above) these states will not be constant and the speed of a shock will change in time.
Shock solution¶
For Burgers' equation, as for the traffic model (or any scalar hyperbolic conservation law with a convex flux function), the solution of the Riemann problem will always consist of one wave, which may be either a shock or a rarefaction. The analysis above already yields the solution to the Riemann problem in the case that the resulting wave is a shock. The entropy condition, as we will explain below, indicates that a shock will form when $f'(q_l) > f'(q_r)$; i.e. when $q_l>q_r$. In this case, the solution is
\begin{align*}
q(x,t) =
\begin{cases}
q_l \quad \text{if} \ \ x/t<s \\
q_r \quad \text{if} \ \ x/t>s.
\end{cases}
\end{align*}
Below, we plot the solution of the Burgers' equation for an initial condition that leads to a shock (i.e. with $q_l>q_r$).
interact(burgers_demos.shock(), t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,description='Show characteristics'));
Rarefaction wave¶
In the previous figure with the hump as the initial condition, we observed that a shock formed on the right side of the hump. However, on the left side, the characteristics spread out and will never cross. This part of the solution is therefore a rarefaction wave. This is the kind of behavior we will observe in the solution of the Riemann problem when $q_l<q_r$.
In the next figure, we consider such a Riemann problem. As time evolves, a rarefaction forms following the advection speed $q$ given by the quasi-linear equation. Therefore, for time $t$, the solution must propagate a distance $x=qt$, so the solution along the rarefaction is then $q = x/t$. As $q_l<q_r$, the smallest and largest displacements are given by $q_l t$ and $q_r t$, respectively.
interact(burgers_demos.rarefaction_figure,
t = FloatSlider(min=0,max=10,value=5));
The full rarefaction solution for the Burgers' Riemann problem is then simply given by \begin{align*} q(x,t) = \begin{cases} q_l, \quad \text{for} \ \ x<q_l t \\ x/t, \quad \text{for} \ \ q_l t \le x \le q_r t \\ q_r, \quad \text{for} \ \ x>q_r t. \end{cases} \end{align*}
As we will see in the next chapter, the rarefaction solution is always a self-similar solution. This means that it can be expressed as a function of the ratio between position and time $q(x,t) = \tilde{q}(x/t)$, so it remains the same when rescaling both $x$ and $t$ by the same factor. This is a consequence of $q$ being the same along the characteristic in the rarefaction fan. In the Burgers' equation, the form of the rarefaction is particularly simple since the advection speed is simply $q$.
Rarefaction solution¶
In the figure below, we plot a solution of the Riemann problem with $q_l<q_r$. In this case, we observe a rarefaction.
interact(burgers_demos.rarefaction(), t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,description='Show characteristics'));
Weak solutions¶
As we mentioned before, the differential form of the equation breaks down in the presence of shocks/discontinuities. However, the integral form of the conservation law remains valid. We will derive a specific form of the integral conservation law that is particularly useful to prove mathematical results. We first integrate the general conservation law $q_t+f(q)_x=0$ from from $x=x_1$ to $x=x_2$ and $t=t_1$ to $t=t_2$
\begin{align} \int_{t_1}^{t_2}\int_{x_1}^{x_2} [q_t+f(q)_x] dx dt = 0, \label{eq:Burgersintclaw} \end{align} This integral can be rewritten in terms of an indicator function $\phi(x,t)$ with compact support in $[x_1,x_2]\times[t_1,t_2]$ (defined to be 1 in this region, zero elsewhere):
\begin{align*} \int_{0}^{\infty}\int_{-\infty}^{\infty} [q_t+f(q)_x]\phi(x,t) dx dt = 0. \end{align*}
We can generalize this result by replacing $\phi(x,t)$ by a smooth function with compact support on some region of the $x-t$ plane. In this case, integrating by parts yield
\begin{align*} \int_{0}^{\infty}\int_{-\infty}^{\infty} [q\phi_t+f(q)\phi_x] dx dt = -\int_{-\infty}^{\infty}q(x,0)\phi(x,0)dx, \end{align*} where now the derivatives are on $\phi(x,t)$ and not on $q$ or $f(q)$, so the equation still makes sense for discontinuous $q$. The function $q(x,t)$ is a weak solution of the conservation law if this equation holds for all functions $\phi(x,t)$ that are continuously differentiable and with compact support (bump functions). Note this rules out the original interval chosen in Eq. \ref{eq:Burgersintclaw} since its $\phi$ would not be smooth. However, we can approximate this function arbitrarily well by a smooth function. Note any weak solution is a solution of any form of the integral conservation law and vice versa.
Weak solutions are thus allowed to include discontinuities. The shock solution presented above, for instance, is a weak solution of the Burgers' equation. After characteristics cross, we cannot have strong solutions of the conservation law, and we must resort to weak solutions. However, one given problem may have multiple weak solutions.
An unphysical weak solution¶
We previously mentioned that for Burgers' equation, when $q_l<q_r$, we expect a rarefaction to be the correct solution. This is because if we smeared out the initial data just slightly then there would be smoothly varying characteristics emerging from each point $x$ and these characteristics would spread out. However, for exactly discontinuous initial data, we could also assume the solution consists of a propagating discontinuity and then use the equation in its conservative form to derive the Rankine-Hugoniot conditions, as done in previous chapters and in (LeVeque, 2002). This is mathematically correct and one obtains the same relation for the shock speed as before, $s = (q_l + q_r)/2$ , so the solution is simply
\begin{align*} q(x,t) = \begin{cases} q_l \quad \text{if} \ \ x/t<s \\ q_r \quad \text{if} \ \ x/t>s. \end{cases} \end{align*}
This unphysical solution is also referred to as an expansion shock, and it is also a weak solution to Burgers' equation. In the next figure, we plot this weak solution.
interact(burgers_demos.unphysical(), t=widgets.FloatSlider(min=0,max=1.0,value=0.5),
which_char=widgets.Checkbox(value=True,description='Show characteristics'));
Note the speed of the characteristics with respect to the shock. The spreading characteristics clearly indicate that the correct solution should be a rarefaction. In order to be able to specify which of the weak solutions is physically correct, we need to derive a mathematical condition from our physical intuition gained from observing the behavior of the characteristics. This condition is referred as the entropy condition.
The entropy condition¶
In the traffic flow chapter, we introduced the entropy condition which determines the correct physical solution among all the possible weak solutions. In the context of the Riemann problem, the entropy condition allows us to determine if the physical solution should involve a shock or a rarefaction. In the case of scalar conservation laws, the entropy condition is quite simple and can be formulated in terms of the flux function of the conservation law as follows: the solution of a scalar Riemann problem will consist of a shock only if
$$ f'(q_l) > f'(q_r). $$ In other words, the solution is a shock only if the characteristics are going into the shock as time progresses. If the characteristics are spreading out/going out of the shock as in the last example, one would naturally expect a rarefaction to be the correct solution. In the case of Burgers' equation, the flux function is $f(q)=q^2/2$, so the correct solution is a shock only if
$$ q_l > q_r, $$ which can be clearly observed in the interactive solution below. The name, "entropy condition", comes from gas dynamics, where the physical entropy must increase in the gas passing through a shock, according to the second law of thermodynamics. A discontinuity that violates this is non-physical. A mathematical version of an "entropy function" can be defined for other conservation laws. More discussion of this and several other formulations and more detailed explanations available in the literature (LeVeque, 2002). In later chapters, we will see how this condition generalizes to systems of conservation laws.
Convexity¶
The flux $f(q) = q^2/2$ for Burgers' equation is a convex function, since $f''(q) = 1 > 0$ for all values of $q$. This means that as $q$ varies between $q_l$ and $q_r$ the the characteristic speed $f'(q) = q$ is either monotonically increasing (if $q_l < q_r$) or monotomically decreasing (if $q_l > q_r$). Because of this the solution is always either a single rarefaction wave or a single shock wave.
The flux $f(\rho) = \rho(1-\rho)$ from the LWR traffic flow model also has the property that $f'(\rho)$ varies monotonically with $\rho$ (in the context of hyperbolic PDEs, a flux function is referred to as convex if either $f''(q)\ge0$ for all $q$ or $f''(q)\le 0$ for all $q$). Since $f''$ is always negative, as we saw in Traffic_flow.html, the correct Riemann solution is a shock if $\rho_l < \rho_r$, or a rarefaction wave if $\rho_l > \rho_r$.
The solution to the Riemann problem can be much more complicated if $f'(q)$ is not monotonically varying between $q_l$ and $q_r$, i.e. if $f''(q)$ changes sign. In this case the Riemann solution can consist of multiple shock and rarefaction waves. The nonconvex case is explored further in Nonconvex_scalar.html.
Interactive solution and examples¶
We can now explore a few more examples that are representative of phenomena we will observe in more complicated systems. Before exploring more complicated examples, we first show a full interactive solution for the Riemann problem for Burgers' equation. The values of the initial conditions and the time can be modified to observe their effect on the characteristic structure and the solution.
burgers.riemann_solution_interact()
burgers_demos.bump_animation(numframes = 50)
Example 3: Three state solutions¶
When three piecewise constant values are given as initial condition, one can solve a Riemann problem across each discontinuity. For instance, consider the Burgers' equation with its initial condition given by
\begin{align*} q_0(x) = \begin{cases} q_l \quad \text{if} \ \ x < -1 \\ q_m \quad \text{if} \ \ -1\le x \le 1 \\ q_r \quad \text{if} \ \ x > 1. \end{cases} \end{align*}
We can decompose this into two Riemann problems, one at $x=-1$ and another one at $x=1$. Each of these two Riemann problems will yield either a shock or a rarefaction. The interesting part is what happens when the generated shocks or rarefactions interact. Following the animation below, let us assume that both Riemann problems produce a right-propagating shock. However, the shock generated at $x=-1$ propagates faster, so it will eventually reach the shock originally generated at $x=1$. At the point in time when the shocks collide, one can simply restate the problem again as a Riemann problem and solve the whole problem analytically. We again use PyClaw to solve the problem and show an animated solution. The first plot shows the solution $q(x)$ evolving through time, while the second one shows the characteristic structure in the $x/t$ plane, where the shocks are shown as wide lines and the characteristics as thin lines. In this figure, one can clearly see the two shocks collide to become a new one. In the second plot, the time is marked by an horizontal dashed line. Note one can also solve this problem analytically. This animation can be viewed in the live notebook on or this webpage.
burgers_demos.triplestate_animation(ql = 4.0, qm = 2.0,
qr = 0.0, numframes = 50)
As one would expect, this can become more complicated when one or both of the waves generated are given by rarefactions. The animation below, or on this webpage, shows how a right-propagating shock can overtake a rarefaction. However, there are consequences for the shock after interacting with the rarefaction; its speed is changed. This is clear in the $x/t$ plane, where the slope of the shock is changed after interacting with the rarefaction. Note the shock is again marked as a wider line.
burgers_demos.triplestate_animation(ql = 4.0, qm = -1.5,
qr = 0.5, numframes = 50)
Analogously, one can also have a rarefaction colliding into a shock. As we can see in the animation below, or this webpage, this can even change the direction in which the shock is propagating. This happens because the velocity of the shock depends on the slopes of the impinging characteristics. When a shock interacts with a rarefaction, the varying slopes of the characteristic curves produce a constant change of the shock propagation velocity. Try modifying the values of $q_l$, $q_m$ and $q_r$ in order to see what other behavior can be observed. Can one observe two interacting rarefactions?
burgers_demos.triplestate_animation(ql = -1, qm = 3.0,
qr = -2, numframes = 50)